function [P_3,T_3,rou_3,Ma_3,a_3]=isentropic_expansion(WAD,kk,P,T,rou,Ma,theta_max,W)
theta=abs(atan(W))+abs(atan(WAD));
Y1=@(x)(theta+((kk+1)/(kk-1))^0.5*atan((((kk-1)/(kk+1))*(Ma^2-1))^0.5)-atan((Ma^2-1)^0.5)-(((kk+1)/(kk-1))^0.5*atan((((kk-1)/(kk+1))*(x^2-1))^0.5)-atan((x^2-1)^0.5)));
d1=1;d2=20;eps=1e-8;
while ((abs(Y1(d1))+abs(Y1(d2)))/2)>=eps
    if theta>theta_max
        d1=20;
        break
    end
    d_t=(d2-d1)/2+d1;
    if (Y1(d1)*Y1(d_t)<0)
        d2=d_t;
    else
        d1=d_t;
    end
    if abs(d1-d2)<1e-12
        'Looping! Help!'
        pause
    end
end
Ma_3=d1;
P_3=((1+((kk-1)/2)*Ma^2)/(1+((kk-1)/2)*Ma_3^2))^(kk/(kk-1))*P;
rou_3=((1+((kk-1)/2)*Ma^2)/(1+((kk-1)/2)*Ma_3^2))^(1/(kk-1))*rou;
T_3=((1+((kk-1)/2)*Ma^2)/(1+((kk-1)/2)*Ma_3^2))*T;
a_3=(kk*P_3/rou_3)^0.5;
end